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Maximum likelihood estimation is applied to the determi- 
nation of an unknown quantum measurement. The measuring 
apparatus performs measurements on many different quan- 
tum states and the positive operator-valued measures govern- 
ing the measurement statistics are then inferred from the col- 
lected data via Maximum-likelihood principle. In contrast to 
the procedures based on linear inversion, our approach always 
provides physically sensible result. We illustrate the method 
on the case of Stern-Gerlach apparatus. 

PACS numbers: 03.65. Bz, 03.67.-a 



I. INTRODUCTION 

Let us imagine that we possess an apparatus which per- 
forms some measurement on certain quantum mechanical 
system such as spin of electron. We do not know which 
measurement is associated with the device and we would 
like to find it out. 

Obviously, the path to follow here is to perform set 
of measurements on various known quantum states and 
then estimate the unknown measurement from the col- 
lected data. Such estimation strategy belongs to the 
broad class of the quantum reconstruction procedures 
which attracted considerable attention recently. The 
quantum state reconstruction has been widely studied 
and now represents a well established technique in many 
branches of quantum physics (for a review, see, e.g. 
The estimation of the quantum mechanical processes de- 
scribing input-output transformations of quantum de- 
vices has been discussed in [|| and the problem of com- 
plete characterization of arbitrary measurement process 
has been recently addressed in W. 

Suppose that the apparatus can respond with k differ- 
ent measurement outcomes. As is well known from the 
theory of quantum measurement JbJ , such device is com- 
pletely characterized by k positive operator valued mea- 
sures (POVM) 11/ which govern the measurement statis- 
tics, 



= i. 



(3) 



pim = Tr[II/p m ], 



(1) 



where g m denotes density matrix of the quantum state 
subject to the measurement, p; m denotes the probability 
that the apparatus would respond with outcome II; to 
the quantum state g m , and Tr stands for the trace. The 
POVMs are positive semi-definite Hermitian operators, 



IT, > 0, 

which decompose the unit operator, 



(2) 



1=1 



The condition (g) ensures that pi m > and (|J) follows 
from the requirement that total probability is normalized 

to unity, Yli=iPlm = l - 

In order to determine the POVMs, one should mea- 
sure on different known quantum states g m and then es- 
timate the POVMs II; from the acquired statistics. Let 
Fi m denote the total number of detections of II/ for the 
measurements performed on the quantum state g m . As- 
suming that the theoretical detection probability pi m (|l|) 
can be replaced with relative frequency, we may write 



N 



Trp/Pm] = II;. 



Fir, 



E*=i F r 



(4) 



where TV is the dimension of the Hilbert space on which 
the operators II/ act. This establishes a system of lin- 
ear equations for the unknown elements of the opera- 
tors Ht,ij, which may easily be solved if sufficient amount 
of data is available. This approach is a direct analogue 
of linear reconstruction algorithms devised for quantum 
state reconstruction. The linear inversion is simple and 
straightforward, but it has also one significant disad- 
vantage. The linear procedure cannot guarantee the re- 
quired properties of II/, namely the conditions (g). Con- 
sequently, the linear estimation may lead to unphysical 
POVMs, predicting negative probabilities pi m for cer- 
tain input quantum states. To avoid such problems, one 
should resort to more sophisticated nonlinear reconstruc- 
tion strategy. 

In this paper we show that the Maximum-likelihood 
(ML) estimation is suitable and can be successfully used 
for the calibration of the measuring apparatus. ML es- 
timation has been recently applied to reconstruction of 
quantum states and quantum processes (complete 
positive maps between density matrices) ||. Here we 
employ it to reconstruct an unknown quantum measure- 
ment thereby demonstrating again the remarkable ver- 
satility and usefulness of ML estimation. The generic 
formalism is developed in Sec. II and an illustration of 
our method is provided in Sec. Ill, where we estimate a 
measurement performed by Stern-Gerlach apparatus. 



II. MAXIMUM-LIKELIHOOD ESTIMATION 

The estimated operators 11/ maximize the likelihood 
functional 



1 



k M 



\ |[ (Tr[n^ m ]) " , 
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(5) 



where 



Jim 



Fl r 



J'=l m'=l 



(6) 



is the relative frequency and M is the number of differ- 
ent quantum states g m used for the reconstruction. The 
maximum of the likelihood functional (j^ has to be found 
in the sub-space of physically allowed operators II; . We 
can decompose each operator 11; as 



JY 



til =^2r tq \cj)i q ){(t)i q \ 

q=l 



(7) 



where ri q > are the eigenvalues of IL and \4>i q ) are 
corresponding orthonormal eigenstates. The maximum 
of £[{11;}] can be found from the extremum conditions. 
It is convenient to work with the logarithm of the origi- 
nal likelihood functional and the constraint (gj) has to be 
incorporated by introducing a Hcrmitian matrix of La- 
grange multipliers A<y = A*^. The extremum conditions 
then read 
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d(<f>i q \ 
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Thus we immediately find 

n q \<t>iq) = Rin q \(f>iq), 

where 

M f 
5 \ — 1 \ Jim ~ 
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and 



N 



(8) 
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(11) 



Let us now multiply (^J) by {4>i q \ from the right and sum 
over q. Thus we obtain 



n, = R,n 



l-u-i- 



(12) 



On averaging ( |l2|) and its Hermitean conjugate counter- 
part, we get 



tli = -(Rifli+fliRj). 



(13) 



The matrix of Lagrange multipliers A should be deter- 
mined from the constraint (^). On summing Eq. (|l"~ 
over I, we find 



k M „ 

y y. —9mUi. 
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(14) 



Eqs. (|T|) and @ can be conveniently solved by means 
of repeated iterations. 

If the linear inversion based on Eqs. (Q) provides phys- 
ically sensible result, then the ML estimate agrees with 
this linear reconstruction. To prove it explicitly, let us 
assume that the set of POVMs IT solves the Eqs. (g). 
Thus we have 



fl 



l in 



X)/'=l fl'm 

On inserting this expression into (p^), we obtain 

m \l' = l ) 



(15) 



(16) 



Here the prime indicates that we should sum only over 
those m with nonzero fi m . However, this restriction may 
be dropped. If /;,„ = pi m = 0, then g m Tli = and the 
addition of zero to the right-hand side of (|l^) changes 
nothing. Thus the set of k equations (|l6|) reduces to 



A I 



lm Qn 



(17) 



m=l 1=1 



which is the formula for the operator of Lagrange multi- 
pliers valid when the measured data are compatible with 
some set of POVMs. Notice that A is positive definite. 

The differences between linear reconstructions and ML 
estimation occur if the experimental data are not compat- 
ible with any physically allowed set of POVMs. The pro- 
cedure of ML estimation may be interpreted as a synthe- 
sis of information from mutually incompatible observa- 
tions. The ML can correctly handle noisy data and pro- 
vides reliable estimates in cases when linear algorithms 
fail. 

Notice that the operators Ri contain the inversion of 
the matrix A. The reconstruction is possible only on such 
subspace of the total Hilbert space where the inversion 
A -1 exists. This restriction can easily be understood if 
we make use of Eq. (p^7|). The experimental data con- 
tain only information on the Hilbert subspace probed by 
the density matrices g m and the reconstruction of the 
POVMs must be restricted to this subspace. 

One could complain that it is not certain that the pos- 
itive definitcness of POVMs 11/ is preserved during iter- 
ations based on Eqs. ( p^ ) and (|l4|). We can, however, 
avoid such complaints by devising an iterative algorithm 
which exactly satisfies the constraints and (||) at each 
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iteration step. We observe that we can formally decom- 
pose the POVMs as 



% = b\b u 



where 



9=1 



(18) 



(19) 



and \q) is some chosen orthonormal basis, (q\q') = 5 qq >. 
From the equations (H) we can derive 



Di = DiR 



(20) 



The constraint (JsJ) provides formula for the operator of 
Lagrange multipliers, 



where G is positive operator 
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(21) 



(22) 



Upon solving ( pl| ) we get A = G 1 / 2 . We fix the branch 
of the square root of G by requiring that A should be 
positive definite operator. We can factorize the matrix 
G as G = U^AU where U is unitary matrix and A is 
diagonal matrix containing eigenvalues of G. We define 



A 1/2 = diag(, 



and we can write 



A = U r k l/2 U. 



(23) 



(24) 



The advantage of the iterative procedure based on (|20| ) 
and (^4|) is that both conditions (||) and (||) are exactly 
fulfilled at each iteration step. The disadvantage of this 
approach is the greater numerical complexity in compar- 
ison to iterations based on (|l2] ) and (fTi]), because we 
must calculate the eigenvalues of the matrix G at each 
iteration step. 

The determination of the quantum measurement can 
simplify considerably if we have some a-priori informa- 
tion about the apparatus. For example, if we know that 
we deal with a photodetector, then we have to estimate 
only a single parameter, the absolute photodetection ef- 
feciency rj |J . Here we briefly consider a broader class of 
phase- insensitive detectors which are sensitive only to the 
number of photons in a single mode of electromagnetic 
field. The POVMs describing phase-insensitive detector 
are all diagonal in the Fock basis, 

IL = },ri n \ri)(n\ 



and the ML estimation reduces to the determination of 
the eigenvalues r;„ > 0. The extremum Eqs. (|l3| ) and 
(|l4|) simplify to 
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(25) 



Instead of solving the extremum equations, one may 
directly search for the maximum of £[{II;}] with the 
help of downhill-simplex algorithm [Q. To implement 
this algorithm successfully, it is necessary to use a min- 
imal parametrization. If we deal with N level system, 
then each 11/ is parametrized by A^ 2 real numbers. Since 
the constraint (||) allows us to determine one POVM in 
terms of remaining fc — 1 ones, the number of independent 
real parameters reads N 2 (k — 1). Furthermore we may 
take advantage of the Cholesky decomposition, 



= cjc h 



(26) 



where Ci is lower triangular matrix with real elements on 
its main diagonal. The parametrization (|26| ) is used for 
the first fc — 1 operators, and the last one is calculated 
form (|), 



IT 



k-l 



cjc h 



(27) 



thus achieving minimal parametrization. For each pa- 
rameter set, where £[{II;}] is evaluated, one has to check 
whether ( |27| ) is positive semi-definite. If this does not 
hold, then one sets £[{II;}] = 0, thus restricting the nu- 
merical search of the maximum to the domain of physi- 
cally allowed operators. This domain is a finite volume 
subspace of a A^ 2 (fc — 1) dimensional space. 



III. STERN-GERLACH APPARATUS 

In this section we illustrate the developed formalism 
by means of numerical simulations for Stern-Gerlach ap- 
paratus measuring a spin-1 particle. We choose the three 
eigenstates of the operator of z-component of the spin as 
the basis states, cr z \s z ) — s z \s z ), s z = —1, 0, 1. In this ba- 
sis, the matrix representation of the spin operators reads 
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FIG. 1: Reconstructed POVMs 11; (solid bars) and the 
exact POVMs (hollow bars). The figure displays elements of 
the vectors fi( defined in Eq. (j2s|) - In the simulation, we 
made M = 30 measurements on each of 12 different quantum 
states, which represents altogether 360 measurements. 

In our numerical simulations, we assume 12 different 
pure quantum states: three eigenstates of a z , \ — l z ), 
|0 Z ), and |l z ), and nine superposition states 2 _1 / 2 (| — 
lz) + e*+' |0,», 2- 1 /2(|o 2 ) + e **i |l z )), and 2-V2(| - 1 2 ) + 
e*|l z )), where j = 1,2,3 and ip x = 0, ip 2 = tt/2, and 
ip3 = ir. The measurement on each state is performed A/" 
times thus the total number of measurements is 127V". We 
have performed Monte Carlo simulations of the measure- 
ment and then we have reconstructed the three POVMs 
III characterizing given apparatus. In the simulations, 
the device measures the projection of the spin along x 
axis, and the true POVMs are projectors 

fti = |ix)<i»|, n 2 = |0x}<0x|, n 3 = l - ix)<-ix|, 

(28) 

&x\sx) — s x \s x ). The elements of the POVMs can be 
conveniently collected in a real vector 

II; = (n^n, IT,, 22, n;,33, ReII U 2, ImlT,^, 

Ren W3 , Im±I U3 , Relics, Imn,, 23 ). (29) 



The estimated POVMs as well as the 'true' POVMs (g|) 
are plotted in Fig. 1. The estimate has been obtained 
by iteratively solving the extremum equations and both 
sets (|H), ((jj) and (§0]), (|J) provide identical results. 
The reconstructed operators are in good agreement with 
the exact ones. Equally important is the fact that the 
estimated operators II; meet the constraints (^) and (^) . 

In summary, we have shown how to reconstruct a 
generic quantum measurement with the use of Maximum- 
likelihood principle. Our method guarantees that the es- 
timated POVMs, which fully describe the measuring ap- 
paratus, meet all the required positivity and completness 
constraints. The numerical feasibility of our technique 
has been illustrated by means of numerical simulations 
for Stern-Gerlach apparatus. 
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